home *** CD-ROM | disk | FTP | other *** search
/ CU Amiga Super CD-ROM 15 / CU Amiga Magazine's Super CD-ROM 15 (1997)(EMAP Images)(GB)[!][issue 1997-10].iso / CUCD / Graphics / Ghostscript / source / gshtscr.c < prev    next >
C/C++ Source or Header  |  1997-05-07  |  17KB  |  523 lines
  1. /* Copyright (C) 1993, 1996, 1997 Aladdin Enterprises.  All rights reserved.
  2.   
  3.   This file is part of Aladdin Ghostscript.
  4.   
  5.   Aladdin Ghostscript is distributed with NO WARRANTY OF ANY KIND.  No author
  6.   or distributor accepts any responsibility for the consequences of using it,
  7.   or for whether it serves any particular purpose or works at all, unless he
  8.   or she says so in writing.  Refer to the Aladdin Ghostscript Free Public
  9.   License (the "License") for full details.
  10.   
  11.   Every copy of Aladdin Ghostscript must include a copy of the License,
  12.   normally in a plain ASCII text file named PUBLIC.  The License grants you
  13.   the right to copy, modify and redistribute Aladdin Ghostscript, but only
  14.   under certain conditions described in the License.  Among other things, the
  15.   License requires that the copyright notice and this notice be preserved on
  16.   all copies.
  17. */
  18.  
  19. /* gshtscr.c */
  20. /* Screen (Type 1) halftone processing for Ghostscript library */
  21. #include "math_.h"
  22. #include "gx.h"
  23. #include "gserrors.h"
  24. #include "gsstruct.h"
  25. #include "gxarith.h"
  26. #include "gzstate.h"
  27. #include "gxdevice.h"            /* for gzht.h */
  28. #include "gzht.h"
  29.  
  30. /* Structure descriptors */
  31. private_st_gs_screen_enum();
  32.  
  33. /* GC procedures */
  34. #define eptr ((gs_screen_enum *)vptr)
  35.  
  36. private ENUM_PTRS_BEGIN(screen_enum_enum_ptrs) {
  37.     if ( index < 2 + st_ht_order_max_ptrs )
  38.       { gs_ptr_type_t ret = (*st_ht_order.enum_ptrs)(&eptr->order, sizeof(eptr->order), index-2, pep);
  39.         if ( ret == 0 )    /* don't stop early */
  40.           ret = ptr_struct_type, *pep = 0;
  41.         return ret;
  42.       }
  43.     return (*st_halftone.enum_ptrs)(&eptr->halftone, sizeof(eptr->halftone), index-(2+st_ht_order_max_ptrs), pep);
  44.     }
  45.     ENUM_PTR(0, gs_screen_enum, pgs);
  46. ENUM_PTRS_END
  47.  
  48. private RELOC_PTRS_BEGIN(screen_enum_reloc_ptrs) {
  49.     RELOC_PTR(gs_screen_enum, pgs);
  50.     (*st_halftone.reloc_ptrs)(&eptr->halftone, sizeof(gs_halftone), gcst);
  51.     (*st_ht_order.reloc_ptrs)(&eptr->order, sizeof(gx_ht_order), gcst);
  52. } RELOC_PTRS_END
  53.  
  54. #undef eptr
  55.  
  56. /* Define the default value of AccurateScreens that affects */
  57. /* setscreen and setcolorscreen. */
  58. private bool screen_accurate_screens = false;
  59.  
  60. /* Default AccurateScreens control */
  61. void
  62. gs_setaccuratescreens(bool accurate)
  63. {    screen_accurate_screens = accurate;
  64. }
  65. bool
  66. gs_currentaccuratescreens(void)
  67. {    return screen_accurate_screens;
  68. }
  69.  
  70. /* Define the MinScreenLevels user parameter similarly. */
  71. private uint screen_min_screen_levels = 1;
  72.  
  73. void
  74. gs_setminscreenlevels(uint levels)
  75. {    screen_min_screen_levels = levels;
  76. }
  77. uint
  78. gs_currentminscreenlevels(void)
  79. {    return screen_min_screen_levels;
  80. }
  81.  
  82. /*
  83.  * The following implementation notes complement the general discussion of
  84.  * halftone tiles found in gxdht.h.
  85.  *
  86.  * Currently we allow R(') > 1 (i.e., multiple basic cells per multi-cell)
  87.  * only if AccurateScreens is true or if B (the number of pixels in a basic
  88.  * cell) < MinScreenLevels; if AccurateScreens is false and B >=
  89.  * MinScreenLevels, multi-cells and basic cells are the same.
  90.  *
  91.  * To find the smallest super-cell for a given multi-cell size, i.e., the
  92.  * smallest (absolute value) coordinates where the corners of multi-cells
  93.  * lie on the coordinate axes, we compute the values of i and j that give
  94.  * the minimum value of W by:
  95.  *    D = gcd(abs(M'), abs(N)), i = M'/D, j = N/D, W = C / D,
  96.  * and similarly
  97.  *    D' = gcd(abs(M), abs(N')), i' = N'/D', j' = M/D', W' = C / D'.
  98.  */
  99.  
  100. /* Compute the derived values of a halftone tile. */
  101. void
  102. gx_compute_cell_values(gx_ht_cell_params_t *phcp)
  103. {    const int M = phcp->M, N = phcp->N, M1 = phcp->M1, N1 = phcp->N1;
  104.     const uint m = any_abs(M), n = any_abs(N);
  105.     const uint m1 = any_abs(M1), n1 = any_abs(N1);
  106.     const ulong C = phcp->C = (ulong)m * m1 + (ulong)n * n1;
  107.     const int D = igcd(m1, n);
  108.     const int D1 = igcd(m, n1);
  109.     const uint W = phcp->W = C / D;
  110.     const uint W1 = phcp->W1 = C / D1;
  111.  
  112.     phcp->D = D, phcp->D1 = D1;
  113.     /* Compute the shift value. */
  114.     /* If M1 or N is zero, the shift is zero. */
  115.     if ( M1 && N )
  116.       {    int h = 0, k = 0, dy = 0;
  117.         int shift;
  118.  
  119.         while ( dy != D )
  120.           if ( dy > D )
  121.             { if ( M1 > 0 ) ++k;
  122.               else --k;
  123.               dy -= m1;
  124.             }
  125.           else
  126.             { if ( N > 0 ) ++h;
  127.               else --h;
  128.               dy += n;
  129.             }
  130.         shift = h * M + k * N1;
  131.         /* We just computed what amounts to a right shift; */
  132.         /* what we want is a left shift. */
  133.         phcp->S = imod(-shift, W);
  134.       }
  135.     else
  136.       phcp->S = 0;
  137.     if_debug12('h', "[h]MNR=(%d,%d)/%d, M'N'R'=(%d,%d)/%d => C=%lu, D=%d, D'=%d, W=%u, W'=%u, S=%d\n",
  138.            M, N, phcp->R, M1, N1, phcp->R1,
  139.            C, D, D1, W, W1, phcp->S);
  140. }
  141.  
  142. /* Forward references */
  143. private int pick_cell_size(P6(gs_screen_halftone *ph,
  144.   const gs_matrix *pmat, ulong max_size, uint min_levels, bool accurate,
  145.   gx_ht_cell_params_t *phcp));
  146.  
  147. /* Allocate a screen enumerator. */
  148. gs_screen_enum *
  149. gs_screen_enum_alloc(gs_memory_t *mem, client_name_t cname)
  150. {    return gs_alloc_struct(mem, gs_screen_enum, &st_gs_screen_enum, cname);
  151. }
  152.  
  153. /* Set up for halftone sampling. */
  154. int
  155. gs_screen_init(gs_screen_enum *penum, gs_state *pgs,
  156.   gs_screen_halftone *phsp)
  157. {    return gs_screen_init_accurate(penum, pgs, phsp,
  158.                        screen_accurate_screens);
  159. }
  160. int
  161. gs_screen_init_memory(gs_screen_enum *penum, gs_state *pgs,
  162.   gs_screen_halftone *phsp, bool accurate, gs_memory_t *mem)
  163. {    int code =
  164.       gs_screen_order_init_memory(&penum->order, pgs, phsp, accurate, mem);
  165.  
  166.     if ( code < 0 )
  167.       return code;
  168.     return
  169.       gs_screen_enum_init_memory(penum, &penum->order, pgs, phsp, mem);
  170. }
  171.  
  172. /* Allocate and initialize a spot screen. */
  173. /* This is the first half of gs_screen_init_accurate. */
  174. int
  175. gs_screen_order_init_memory(gx_ht_order *porder, const gs_state *pgs,
  176.   gs_screen_halftone *phsp, bool accurate, gs_memory_t *mem)
  177. {    gs_matrix imat;
  178.     ulong max_size = pgs->ht_cache->bits_size;
  179.     int code;
  180.  
  181.     if ( phsp->frequency < 0.1 )
  182.       return_error(gs_error_rangecheck);
  183.     gs_deviceinitialmatrix(gs_currentdevice(pgs), &imat);
  184.     code = pick_cell_size(phsp, &imat, max_size,
  185.                   screen_min_screen_levels, accurate,
  186.                   &porder->params);
  187.     if ( code < 0 )
  188.       return code;
  189.     gx_compute_cell_values(&porder->params);
  190.     /* NOTE: patching the next line to 'false' forces all screens */
  191.     /* to be strip halftones, even if they are small. */
  192.     if ( porder->params.W1 <= max_size / bitmap_raster(porder->params.W) )
  193.       { /*
  194.          * Allocate an order for the entire tile, but only sample one
  195.          * strip.  Note that this causes the order parameters to be
  196.          * self-inconsistent until gx_ht_construct_spot_order fixes them
  197.          * up: see gxdht.h for more information.
  198.          */
  199.         code = gx_ht_alloc_order(porder, porder->params.W,
  200.                      porder->params.W1, 0,
  201.                      porder->params.W * porder->params.D,
  202.                      mem);
  203.         porder->height = porder->orig_height = porder->params.D;
  204.         porder->shift = porder->orig_shift = porder->params.S;
  205.       }
  206.     else
  207.       { /* Just allocate the order for a single strip. */
  208.         code = gx_ht_alloc_order(porder, porder->params.W,
  209.                      porder->params.D, porder->params.S,
  210.                      porder->params.W * porder->params.D,
  211.                      mem);
  212.       }
  213.     if ( code < 0 )
  214.       return code;
  215.     return 0;
  216. }
  217.  
  218. /*
  219.  * Given a desired frequency, angle, and minimum number of levels, a maximum
  220.  * cell size, and an AccurateScreens flag, pick values for M('), N('), and
  221.  * R(').  We want to get a good fit to the requested frequency and angle,
  222.  * provide at least the requested minimum number of levels, and keep
  223.  * rendering as fast as possible; trading these criteria off against each
  224.  * other is what makes the code complicated.
  225.  *
  226.  * We compute trial values u and v from the original values of F and A.
  227.  * Normally these will not be integers.  We then examine the 4 pairs of
  228.  * integers obtained by rounding each of u and v independently up or down,
  229.  * and pick the pair U, V that yields the closest match to the requested
  230.  * F and A values and doesn't require more than max_size storage for a
  231.  * single tile.  If no pair
  232.  * yields an acceptably small W, we divide both u and v by 2 and try again.
  233.  * Then we run the equations backward to obtain the actual F and A.
  234.  * This is fairly easy given that we require either xx = yy = 0 or
  235.  * xy = yx = 0.  In the former case, we have
  236.  *    U = (72 / F * xx) * cos(A);
  237.  *    V = (72 / F * yy) * sin(A);
  238.  * from which immediately
  239.  *    A = arctan((V / yy) / (U / xx)),
  240.  * or equivalently
  241.  *    A = arctan((V * xx) / (U * yy)).
  242.  * We can then obtain F as
  243.  *    F = (72 * xx / U) * cos(A),
  244.  * or equivalently
  245.  *    F = (72 * yy / V) * sin(A).
  246.  * For landscape devices, we replace xx by yx, yy by xy, and interchange
  247.  * sin and cos, resulting in
  248.  *    A = arctan((U * xy) / (V * yx))
  249.  * and
  250.  *    F = (72 * yx / U) * sin(A)
  251.  * or
  252.  *    F = (72 * xy / V) * cos(A).
  253.  */
  254. /* ph->frequency and ph->angle are input parameters; */
  255. /* the routine sets ph->actual_frequency and ph->actual_angle. */
  256. private int
  257. pick_cell_size(gs_screen_halftone *ph, const gs_matrix *pmat, ulong max_size,
  258.   uint min_levels, bool accurate, gx_ht_cell_params_t *phcp)
  259. {    const bool landscape = (pmat->xy != 0.0 || pmat->yx != 0.0);
  260.     /* Account for a possibly reflected coordinate system. */
  261.     /* See gxstroke.c for the algorithm. */
  262.     const bool reflected = pmat->xy * pmat->yx > pmat->xx * pmat->yy;
  263.     const int reflection = (reflected ? -1 : 1);
  264.     const int rotation =
  265.       (landscape ? (pmat->yx < 0 ? 90 : -90) : pmat->xx < 0 ? 180 : 0);
  266.     const double f0 = ph->frequency, a0 = ph->angle;
  267.     const double T =
  268.       fabs((landscape ? pmat->yx / pmat->xy : pmat->xx / pmat->yy));
  269.     gs_point uv0;
  270. #define u0 uv0.x
  271. #define v0 uv0.y
  272.     int rt = 1;
  273.     double f = 0, a = 0;
  274.     double f_best = max_int, a_best = 360, e_best = 1000;
  275.     long w_size = max_size;
  276.     bool better;
  277.  
  278.     /*
  279.      * We need to find a vector in device space whose length is
  280.      * 1 inch / ph->frequency and whose angle is ph->angle.
  281.      * Because device pixels may not be square, we can't simply
  282.      * map the length to device space and then rotate it;
  283.      * instead, since we know that user space is uniform in X and Y,
  284.      * we calculate the correct angle in user space before rotation.
  285.      */
  286.  
  287.     /* Compute trial values of u and v. */
  288.  
  289.     { gs_matrix rmat;
  290.       gs_make_rotation(a0 * reflection + rotation, &rmat);
  291.       gs_distance_transform(72.0 / f0, 0.0, &rmat, &uv0);
  292.       gs_distance_transform(u0, v0, pmat, &uv0);
  293.       if_debug10('h', "[h]Requested: f=%g a=%g mat=[%g %g %g %g] max_size=%lu min_levels=%u =>\n     u=%g v=%g\n",
  294.              ph->frequency, ph->angle,
  295.              pmat->xx, pmat->xy, pmat->yx, pmat->yy,
  296.              max_size, min_levels, u0, v0);
  297.     }
  298.  
  299.     /* Adjust u and v to reasonable values. */
  300.  
  301.     if ( u0 == 0 && v0 == 0 )
  302.       return_error(gs_error_rangecheck);
  303.     while ( (fabs(u0) + fabs(v0)) * rt < 4 )
  304.       ++rt;
  305. try_size:
  306.     better = false;
  307.     { int m0 = (int)floor(u0 * rt + 0.0001);
  308.       int n0 = (int)floor(v0 * rt + 0.0001);
  309.       gx_ht_cell_params_t p;
  310.  
  311.       p.R = p.R1 = rt;
  312.       for ( p.M = m0 + 1; p.M >= m0; p.M-- )
  313.         for ( p.N = n0 + 1; p.N >= n0; p.N-- )
  314.           {    long raster, wt, wt_size;
  315.         double fr, ar, ft, at, f_diff, a_diff, f_err, a_err;
  316.  
  317.         p.M1 = (int)(p.M / T + 0.5);
  318.         p.N1 = (int)(p.N * T + 0.5);
  319.         gx_compute_cell_values(&p);
  320.         if_debug3('h', "[h]trying m=%d, n=%d, r=%d\n", p.M, p.N, rt);
  321.         wt = p.W;
  322.         if ( wt >= max_short )
  323.           continue;
  324.         /* Check the strip size, not the full tile size, */
  325.         /* against max_size. */
  326.         raster = bitmap_raster(wt);
  327.         if ( raster > max_size / p.D || raster > max_long / wt )
  328.           continue;
  329.         wt_size = raster * wt;
  330.  
  331.         /* Compute the corresponding values of F and A. */
  332.  
  333.         if ( landscape )
  334.           ar = atan2(p.M * pmat->xy, p.N * pmat->yx),
  335.           fr = 72.0 * (p.M == 0 ? pmat->xy / p.N * cos(ar) :
  336.                    pmat->yx / p.M * sin(ar));
  337.         else
  338.           ar = atan2(p.N * pmat->xx, p.M * pmat->yy),
  339.           fr = 72.0 * (p.M == 0 ? pmat->yy / p.N * sin(ar) :
  340.                    pmat->xx / p.M * cos(ar));
  341.         ft = fabs(fr) * rt;
  342.         /* Normalize the angle to the requested quadrant. */
  343.         at = (ar * radians_to_degrees - rotation) * reflection;
  344.         at -= floor(at / 180.0) * 180.0;
  345.         at += floor(a0 / 180.0) * 180.0;
  346.         f_diff = fabs(ft - f0);
  347.         a_diff = fabs(at - a0);
  348.         f_err = f_diff / fabs(f0);
  349.         /*
  350.          * We used to compute the percentage difference here:
  351.          *    a_err = (a0 == 0 ? a_diff : a_diff / fabs(a0));
  352.          * but using the angle difference makes more sense:
  353.          */
  354.         a_err = a_diff;
  355.  
  356.         if_debug5('h', " ==> d=%d, wt=%ld, wt_size=%ld, f=%g, a=%g\n",
  357.               p.D, wt, bitmap_raster(wt) * wt, ft, at);
  358.  
  359.         /*
  360.          * Minimize angle and frequency error within the
  361.          * permitted maximum super-cell size.
  362.          */
  363.  
  364.         { double err = f_err * a_err;
  365.           if ( err > e_best )
  366.             continue;
  367.           e_best = err;
  368.         }
  369.         *phcp = p;
  370.         f = ft, a = at;
  371.         w_size = wt_size, f_best = f_diff, a_best = a_diff;
  372.         better = true;
  373.         if_debug3('h', "*** w_size=%ld, f_best=%g, a_best=%g\n",
  374.               w_size, f_best, a_best);
  375.         if ( f_err <= 0.01 && a_err <= 0.01 )
  376.           goto done;
  377.           }
  378.     }
  379.     if ( phcp->C < min_levels )
  380.       {    /* We don't have enough levels yet.  Keep going. */
  381.         ++rt;
  382.         goto try_size;
  383.       }
  384.     if ( better )
  385.       {    /* If we want accurate screens, continue till we fail. */
  386.         if ( accurate )
  387.           { ++rt;
  388.             goto try_size;
  389.           }
  390.       }
  391.     else
  392.       {    /*
  393.          * We couldn't find an acceptable M and N.  If R > 1,
  394.          * take what we've got; if R = 1, give up.
  395.          */
  396.         if ( rt == 1 )
  397.           return_error(gs_error_rangecheck);
  398.       }
  399.  
  400.     /* Deliver the results. */
  401. done:
  402.     if_debug5('h', "[h]Chosen: f=%g a=%g M=%d N=%d R=%d\n",
  403.           f, a, phcp->M, phcp->N, phcp->R);
  404.     ph->actual_frequency = f;
  405.     ph->actual_angle = a;
  406.     return 0;
  407. #undef u0
  408. #undef v0
  409. }
  410.  
  411. /* Prepare to sample a spot screen. */
  412. /* This is the second half of gs_screen_init_accurate. */
  413. int
  414. gs_screen_enum_init_memory(gs_screen_enum *penum, const gx_ht_order *porder,
  415.   gs_state *pgs, gs_screen_halftone *phsp, gs_memory_t *mem)
  416. {    penum->pgs = pgs;        /* ensure clean for GC */
  417.     penum->order = *porder;
  418.     penum->halftone.rc.memory = mem;
  419.     penum->halftone.type = ht_type_screen;
  420.     penum->halftone.params.screen = *phsp;
  421.     penum->x = penum->y = 0;
  422.     penum->strip = porder->num_levels / porder->width;
  423.     penum->shift = porder->shift;
  424.     /*
  425.      * We want a transformation matrix that maps the parallelogram
  426.      * (0,0), (U,V), (U-V',V+U'), (-V',U') to the square (+/-1, +/-1).
  427.      * If the coefficients are [a b c d e f] and we let
  428.      *    u = U = M/R, v = V = N/R,
  429.      *    r = -V' = -N'/R', s = U' = M'/R',
  430.      * then we just need to solve the equations:
  431.      *    a*0 + c*0 + e = -1    b*0 + d*0 + f = -1
  432.      *    a*u + c*v + e = 1    b*u + d*v + f = 1
  433.      *    a*r + c*s + e = -1    b*r + d*s + f = 1
  434.      * This has the following solution:
  435.      *    Q = 2 / (M*M' + N*N')
  436.      *    a = Q * R * M'
  437.      *    b = -Q * R' * N
  438.      *    c = Q * R * N'
  439.      *    d = Q * R' * M
  440.      *    e = -1
  441.      *    f = -1
  442.      */
  443.     { const int M = porder->params.M, N = porder->params.N,
  444.         R = porder->params.R;
  445.       const int M1 = porder->params.M1, N1 = porder->params.N1,
  446.         R1 = porder->params.R1;
  447.       double Q = 2.0 / ((long)M * M1 + (long)N * N1);
  448.       penum->mat.xx = Q * (R * M1);
  449.       penum->mat.xy = Q * (-R1 * N);
  450.       penum->mat.yx = Q * (R * N1);
  451.       penum->mat.yy = Q * (R1 * M);
  452.       penum->mat.tx = -1.0;
  453.       penum->mat.ty = -1.0;
  454.     }
  455.     if_debug7('h', "[h]Screen: (%dx%d)/%d [%f %f %f %f]\n",
  456.           porder->width, porder->height, porder->params.R,
  457.           penum->mat.xx, penum->mat.xy,
  458.           penum->mat.yx, penum->mat.yy);
  459.     return 0;
  460. }
  461.  
  462. /* Report current point for sampling */
  463. int
  464. gs_screen_currentpoint(gs_screen_enum *penum, gs_point *ppt)
  465. {    gs_point pt;
  466.     int code;
  467.     if ( penum->y >= penum->strip )        /* all done */
  468.     {    gx_ht_construct_spot_order(&penum->order);
  469.         return 1;
  470.     }
  471.     /* We displace the sampled coordinates very slightly */
  472.     /* in order to reduce the likely number of points */
  473.     /* for which the spot function returns the same value. */
  474.     if ( (code = gs_point_transform(penum->x + 0.501, penum->y + 0.498, &penum->mat, &pt)) < 0 )
  475.         return code;
  476.     if ( pt.x < -1.0 )
  477.         pt.x += ((int)(-ceil(pt.x)) + 1) & ~1;
  478.     else if ( pt.x >= 1.0 )
  479.         pt.x -= ((int)pt.x + 1) & ~1;
  480.     if ( pt.y < -1.0 )
  481.         pt.y += ((int)(-ceil(pt.y)) + 1) & ~1;
  482.     else if ( pt.y >= 1.0 )
  483.         pt.y -= ((int)pt.y + 1) & ~1;
  484.     *ppt = pt;
  485.     return 0;
  486. }
  487.  
  488. /* Record next halftone sample */
  489. int
  490. gs_screen_next(gs_screen_enum *penum, floatp value)
  491. {    ht_sample_t sample;
  492.     int width = penum->order.width;
  493.     if ( value < -1.0 || value > 1.0 )
  494.         return_error(gs_error_rangecheck);
  495.     /* The following statement was split into two */
  496.     /* to work around a bug in the Siemens C compiler. */
  497.     sample = (ht_sample_t)(value * max_ht_sample);
  498.     sample += max_ht_sample;    /* convert from signed to biased */
  499. #ifdef DEBUG
  500. if ( gs_debug_c('H') )
  501.    {    gs_point pt;
  502.     gs_screen_currentpoint(penum, &pt);
  503.     dprintf6("[H]sample x=%d y=%d (%f,%f): %f -> %u\n",
  504.          penum->x, penum->y, pt.x, pt.y, value, sample);
  505.    }
  506. #endif
  507.     penum->order.bits[penum->y * width + penum->x].mask = sample;
  508.     if ( ++(penum->x) >= width )
  509.         penum->x = 0, ++(penum->y);
  510.     return 0;
  511. }
  512.  
  513. /* Install a fully constructed screen in the gstate. */
  514. int
  515. gs_screen_install(gs_screen_enum *penum)
  516. {    gx_device_halftone dev_ht;
  517.  
  518.     dev_ht.rc.memory = penum->halftone.rc.memory;
  519.     dev_ht.order = penum->order;
  520.     dev_ht.components = 0;
  521.     return gx_ht_install(penum->pgs, &penum->halftone, &dev_ht);
  522. }
  523.